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The density and correlations of topological defects are investigated numerically in a model of a 
d — 2 elastic medium subject to a periodic quenched random potential. The computed density of 
defects decreases approximately exponentially with the defect core energy. Comparing the defect- 
free ground state with the ground state with defects, it is found that the difference is described by 
string-like excitations, bounded by defect pairs, which have a fractal dimension of 1.250 ± 0.003. At 
zero temperature, the disorder-induced defects screen the interaction of introduced vortex pairs. 
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A diverse set of physical phenomena, including vortex 
lattices in superconductors, incommensurate charge den- 
sity waves, and crystal growth on a disordered substrate, 
have been modeled as elastic media subject to a pinning 
potential due to quenched disorder . A possible dif- 
ficulty in applying these purely elastic models is their 
failure to take into account account the possible effect 
of defects in the elastic medium |QM§. In general, the 
importance of defects have been difficult to analyze an- 
alytically, though recent numerical simulations and ana- 
lytical arguments have indicated that defects significantly 
change the d = 2 system, but have less effect in d = 3 
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In this paper, I present results on the effects of scalar 
defects on the ground state of a model of a d — 2 elastic 
medium subject to a disordered pinning potential. Us- 
ing variants of combinatorial techniques that have been 
applied |^-||] to the study of the elastic model without 
defects, it is possible to exactly calculate the finite-size 
ground states of models that allow for defects in the 
scalar displacement variable. By comparing configura- 
tions with and without defects, the change in ground 
state due to the introduction of defects is found to be 
confined to "strings" that connect defects. The fractal 
dimension of these strings is computed to be 1.250(3). 
The effects of the defects on the long-range response are 
calculated by the introduction of a fixed pair of defects 
into the elastic medium with defects. As the system size 
increases, the cost of introducing the fixed pair goes to a 
constant (rather than growing logarithmically), indicat- 
ing that the defects screen the long-range interaction at 
long scales. 

In models of d-dimensional elastic media with scalar 
displacement subject to a pinning potential periodic in 
the direction of displacement, there arc two length scales 
that separate distinct behaviors WM. For lengths less 
than the pinning length £p, the elastic energies for defor- 
mation are typically larger than the pinning energies, so 
that the displacement in a region of volume £p can be rep- 
resented as a single scalar variable h, with the pinning en- 
ergy periodic in h. For length scales greater than £p, the 



pinning energy dominates and there are many metastable 
states. In a medium that allows for the creation of de- 
fects (in the d = 2 case considered in this paper, these 
defects are point vortices), there is a length scale £y de- 
scribing the typical separation of vortices. The vortices 
themselves are described by their location to within £p 
GIB. In this work, the lattice constant corresponds to 
the length scale £p and the distance £y is controlled by 
a vortex core energy E c . 

A model for interfaces and defects in d = 2 can be 
based upon the properties of matchings and their height 
representations. Given an undirected graph G = (V, D), 
with vertices V and edges D, a matching M C D is a 
set of edges such that each vertex belongs to at most one 
edge in M. The graphs of interest here will be bipartite 
(with subsets of vertices A and B and all edges connect- 
ing A- and B-vertices) , with equal numbers of vertices in 
each sublattice, and allow for complete matchings, where 
each vertex is a member of exactly one edge. In partic- 
ular, consider a hexagonal lattice with L x L two-vertex 
unit cells. There exists a one-to-one mapping between 
complete matchings on this lattice, minimally frustrated 
states of an antiferromagnetic Ising model on a triangu- 
lar lattice, and a solid-on-solid height representation ||. 
Fig. [l] indicates the mapping between matchings and a 
height representation. In the case of a partial matching, 
an equal number of A- and i?-vertices are unmatched. 
The uncovered vertices correspond to defects (vortices 
or screw dislocations) in the height representation. The 
rules for calculating the height can still be applied locally, 
but are inconsistent on any loop containing an unequal 
number of A and B defects. Unmatched vertices have a 
positive (negative) sign, if the vertex is a member of the 
B (A) sublattice, respectively, as indicated in Fig. |l[ 

Consider the case of a complete matching M, where 
there is a well-defined interface corresponding to M. In 
the absence of quenched disorder, the energy of the in- 
terface is independent of the height configuration, so all 
matchings appear with equal weight, and the thermally- 
averaged height-height correlations < [h(r) — h(0)} 2 > 
grows as ~ ln(r) at all temperatures T H. In the 
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presence of quenched disorder, this is no longer true at 
low temperatures. Theoretical calculations |l],[l(| predict 
that the 2+1-dimensional model at low temperature T 
has a correlation function that grows as ~ ln 2 (r). This 
has been confirmed in the T — > limit by numerical 
calculations using the complete matching representation 
and other methods ||,^] . The energy can be calculated in 
the matching representation as E = — YleEM w e> wnere 
w e are randomly chosen weights associated with each 
edge (dimer.) Maximizing the total weight minimizes the 
energy. When the corresponding interface in the height 
representation is translated by 3 units in the height direc- 
tion, the dimer configuration is unchanged, so that the 
energy is periodic in global height shifts. 




FIG. 1. Rules for assigning heights to matchings on a 
hexagonal lattice. A- (B-) vertices are shown as filled (un- 
filled) circles and matched edges are heavy line segments. All 
but the two circled vertices are matched. Heights h are de- 
fined at the centers of the hexagons. When crossing a matched 
edge from an A to a B vertex in the ccw direction about the 
A vertex, h increases by 2; if the edge is unmatched, h de- 
creases by 1. A locally consistent definition of h is shown on 
the left side of the figure. The right side of the figure ex- 
hibits defects. The unmatched vertex on the A (B) sublattice 
correspond to a negative (positive) vortex. The dashed line 
indicates a possible path for a defect string that connects two 
defects; if this were the defect string, then in the ground state 
with no defects, the unmatched edges along this path would 
become matched and the matched edges would be removed. 

In the case of an incomplete matching, there is no 
longer a uniquely defined height variable corresponding 
to displacements of the d — 2 elastic medium. An energy 
for the dimer configuration can then be assigned which is 
consistent with (1) a local pinning energy where vertices 
are covered by dimers in M and (2) a defect energy for 
vertices that are not the endpoints of dimers in M: 



E = N r E r 



(1) 



with defects that have an associated core energy. Over 
local regions where the height is well defined, the energy 
is still periodic in the height (displacement) variable. 

Ground-state configurations were found using two al- 
gorithms that solve maximum- weight bipartite matching 
(MWBM.) The first algorithm uses a heuristic developed 
by C. Zeng [ Il[, implemented in the LEDA library of al- 
gorithms |l2 1, to directly determine the (partial) match- 
ing that maximizes the sum in Eqn. ([!]). The second 
algorithm, using almost identical processor time and less 
memory, is based upon finding a maximum-weight com- 
plete matching in a modified graph, where vertices are 
duplicated and extra edges of weight 2E C and of weight 
zero are added to the original edge set. A complete 
matching algorithm such as the cost-scaling assignment 
(CSA) algorithm by Kennedy and Goldberg |l3| can then 
be used to find a complete matching on this augmented 
graph. Matched vertices in the augmented graph which 
correspond to edges in the original graph give the par- 
tial matching in the original graph that minimizes E. 
These algorithms take 340 s on a 500 MHz DEC Alpha 
to find the ground state of a system with 294 912 vertices 
(384 x 384 unit cells). 

In addition, one can study the result of introducing 
a single defect pair, at given locations, into a complete 
matching M (a defect- free medium.) This is done by 
uncovering an A- and a B-vertex and arranging \M\ — 1 
dimers to cover the remaining vertices. This arrange- 
ment can be found by solving a shortest paths problem 
p in a directed graph G' which is determined by M and 
G. An undirected edge e S G is replaced with a directed 
edge e' G G' from the A- to a B-vertex if e S M, other- 
wise e' is directed from the B- to the A-node. A directed 
path p following the edges in G' is assigned the energy 
change AE(p) = J2 



where N c is the number of non-covered vertices and E c 
is the core cost of a defect. Minimizing E in Eqn. (|l|) 
gives the T = configuration for the elastic medium 



where w' e — w e if e € M and 
w e = —w e for e ^ M. Minimal energy change paths 
starting from an introduced positive defect to any other 
vertex then give the minimal energy excitation due to a 
defect pair introduced at the endpoints of the path. Such 
paths can be determined by a shortest path (SP) algo- 
rithm that allows for negative weights; for this purpose, 
the Goldberg-Razik algorithm was used @. Note that 
this algorithm is distinct from the shortest paths algo- 
rithms used to study the directed polymer problem [fl5| , 
where the edge weights can all be made positive by a 
uniform shift without affecting the paths. 

Finally, the combination of the introduction of a fixed 
pair and defects with a chosen core energy was studied 
using a MWBM algorithm on a graph G" = (V, D U {z}) 
with z an external edge connecting an A-node and a B- 
node at a separation r^e/ = (L/2, L/2) (in the lattice unit 
vector representation), with w z — oo, so that z always 
introduces a pair of defects separated as far as possible 
in a finite system. The MWBM algorithm then gives the 
minimum E for the introduced defects, in the presence 
of disorder- induced defects. 
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Simulations using MWBM without an introduced de- 
fect were performed for a variety of system sizes and de- 
fect core energies. The edge weights (pinning energies) 
were chosen from a uniform distribution in the range 
[0,1), with a discrete resolution of 10 -4 . Defect den- 
sities were computed by comparing \M\ with the number 
of edges in a complete matching (a decrease in dirtier 
number of one gives a defect pair.) Complete matching 
(E c = oo) configurations were compared with the partial 
matchings for identical disorder to determine the changes 
due to defects. For the largest system size (384 x 384 unit 
cells), 400 samples were studied. Simulations for intro- 
duced defect pairs using the SP algorithm were performed 
for systems up to 768 x 768 unit cells (> 4 360 samples.) 

In the MWBM calculation with finite core energy, de- 
fects are readily identified as unmatched vertices (+ if 
a B-vertex, — if an A- vertex.) The total defect density 
p = (N + + NS)/2WL is fit to within statistical errors by 
a simple exponential, p tx cxp(— Ec/Eq), with Eq = 0.45, 
for p^ 1 !" 1 > 10. This form is consistent with the mini- 
mum core energy E™ 1U (L) that typically excludes intro- 
duced defects on a length scale of L; numerical calcula- 
tion of this quantity using the SP algorithm on complete 
matchings gives Ef ia (L) w (const.) + (0.36) ln(L). 

The ground states with and without defects can be 
directly compared. There are two types of changes ap- 
parent: the introduction of vortices and changes in the 
matching in the regions between vortices. The changes 
are realized as "strings" connecting pairs of defects, as 
shown in Fig. ||. This result is consistent with the work 
of Gingras and Huse Q who argued that in an XY model, 
the phase difference caused by a defect will be confined to 
lines connecting vortices of opposite sign. This confine- 
ment results from strong pinning on scales larger than 
£p; the height difference of ±3 found on a loop enclosing 
a vortex is not spread out uniformly, but takes place over 
the scale £p. 




FIG. 2. The symmetric difference between matchings for 
a ground state with no defects and a ground state for the 
same disorder realization, with a defect energy of E c = 1.2 
(384 x 384 unit cells.) Dimers are included if they belong to a 
matching in one of the ground states, but not both. Vortices 
are at the ends of the defect lines. The defect lines them- 
selves show where the phase change due to the introduction 
of defects is localized. 
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FIG. 3. (a) A plot of the mean end-to-end displacement 
r(E c ) of defect strings vs. the mean length s(E c ) of a de- 
fect path, for core energies E c — 0, 0.2, 0.4,. . . , 1.4, found by 
examining defect strings as shown in Fig. |2j Statistical and 
systematic (from the extrapolation L — > oo) uncertainties are 
less than 1% in each quantity. The straight line indicates a 
fit to the power law r ~ s 1/d f , for df = 1.261(16). (b) A 
finite-size scaling plot of r df /s as a function of s/lA from 
spanning trees generated by the SP algorithm on complete 
matchings. In this case, precise collapse over two decades 
with a linear abscissa gives the fitted value of df as 1.250(3). 

The strings corresponding to the difference between 
ground states with and without defects were character- 
ized by the string length s (number of edges in the string) 
and the end-to-end distance r. The mean of the end-to- 
end displacement distance of the paths, r(E c ), and the 
mean length of the paths, s(E c ), depends on the core 
energy. The dimension df was computed assuming 

r(E c )~ a V*t(E e ). (2) 

The quantities r(L, E c ) and s(L, E c ) were simply com- 
puted by averaging over all defect paths at fixed L and 
E c . The limits s(L — » oo,E c ) = s(E c ) and r(L — » 
oo,E c ) — r(E c ) were found by scaling fits. For example, 
a plot of r(L, E c )/r(oo, E c ) vs. L/r(oo, E c ) was made for 
trial values of r(oo, E c ) until the data collapsed to a single 
line. The values of s(L, E c ) were scaled using a presumed 
scaling form s(L,E c )/s(oc,E c ) ~ L/s(oo, E^^f , for 
trial values of s(oo, E c ) and df. The final value of df was 
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then most accurately read from a plot of r(oo,E c ) vs. 
s(oo,E c ), as shown in Fig. 3(a). The single defect paths 
caused by the two vortices placed in a complete match- 
ing, computed using the CSA and SP algorithms, were 
also studied. In a single sample, s can be calculated for 
all possible values of r; the finite-size scaling is shown in 
Fig. 3(b). The fractal dimensions for the multiple defect 
and single introduced defect, df = 1.261(16), 1.250(3), 
respectively, are numerically consistent. 

The SP algorithm gives the minimum energy excita- 
tions for a pair of defects at all separations, in the ab- 
sence of other defects. Applying MWBM to the ex- 
tended graph G" allows one to study a pair of defects 
at fixed separation in the presence of defects. The cost 
AEdef to introduce a defect pair is given by the dif- 
ference between the ground state energies for G and 
G" at the same core energy. Introducing a pair of de- 
fects, if other defects can be ignored, costs an energy 
2E C (for the defect cores) and a separation energy, log- 
arithmic in the separation so that at small values 
of x — p x l 2 L, AE de f = 2E C + c + 27T7p s ln(L), where 
p s is the long-distance stiffness associated with a single 
core, c is a constant, and 7 is a constant given by the 
sample geometry. For large values of x, the introduced 
pair cost approaches a constant, if the disorder- induced 
defects screeen the interaction between the introduced 
defects, with the crossover occuring for p 1 ' 2 L ~ 1 (giv- 
ing AE def = 2E c + c' + 27T7 / 9 s ln( / 9~ 1 / 2 ) for large x.) The 
results of the simulation, as plotted in Fig. are con- 
sistent with this form of defect screening. 
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FIG. 4. A plot of the energy cost AE C to introduce two 
maximally separated defects into a system of size LxL, scaled 
to indicate consistency with the limits described in the text. 
The symbols indicate the core energy E c ; distinct points with 
the same symbol are found by varying L (L = 6, . . . , 192.) Er- 
ror bars are statistical (lcr.) The coefficient yp a was chosen 
for the best scaling collapse. For fixed E c , the defect cost ap- 
proaches a constant in large systems, indicating the screening 
of the defect-defect interaction. 

By applying optimization algorithms to many large 
samples, the lack of long-range elastic behavior in a disor- 
dered medium at T = with defects is confirmed and the 
string dimension df is precisely determined. Note that 



df appears to be distinct from the defect-step dimension 
1.35 ± 0.02 previously reported this difference is sur- 
prising and further confirmation would be useful. The 
dimension df is very close to that of 5/4 for loop-erased 
random walks (LERW's) |16|, and distinct from the di- 
mension (1.22(1) 0; 1.222(3) Jl8|) of minimal-spanning 
trees (if the weights were all positive in G", the search for 
defect strings would give a minimal spanning tree; note 
that the SP algorithm has a loop-cancellation similar to 
that of the LERW.) 
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